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Abstract 

Results  are  presented  of  studies  of  different  types  of  optimal  and  quasi-optimal  Kalman  filters 
based  on  crystal  and  rubidium  oscillators  using  reference  timing  signals  from  the  Motorola  GPS  UT  + 
Oncore  Timing  receiver.  Filter  equations  are  considered  for  different  definitions  of  their  coefficients, 
and  the  filter  output  signal  and  its  statistics  are  investigated  in  time  under  real  conditions.  Various 
Kalman  algorithms  and  corresponding  optimal  filter  structures  intended  for  crystal  and  rubidium 
oscillators  are  discussed.  Three-dimensional  Kalman  filters  intended  for  optimal  estimation  of  time 
error,  frequency  offset,  and  frequency  aging  are  considered  based  on  the  oscillator  signal  model.  One 
application  is  the  synchronization  needs  of  digital  communication  networks  and  metrology.  Results 
are  given  for  the  measured  data  and  estimates  are  compared  with  respect  to  a  quartz  crystal  oscillator. 
Practical  results  are  considered  of  the  Kalman  filter’s  use  in  application  to  an  oven-controlled  quartz 
crystal  oscillator  (OCXO)  with  an  AT-cut  resonator.  Plots  of  the  original  and  filtered  processes  are 
discussed  for  the  different  approaches  to  definition  of  the  Kalman  filter  coefficients.  Estimates  are 
also  given  for  the  filtering  errors  and  the  processing  rate. 


INTRODUCTION 

Fast  and  accurate  optimal  Kalman  filtering  of  the  time  error,  frequency,  and  frequency  offset  of  slaved 
sources  (crystal  and  rubidium)  is  extremely  important  to  create  frequency  and  time  standards  looked  after 
GPS  timing  signals.  It  is  known  that  excellent  accuracy  less  then  10‘12  is  obtained  in  practice 
approximately  through  24  hours,  based  on  a  smoothing  filter.  It  is  the  reason  why  efforts  are  now 
underway  in  direction  of  the  fast  Kalman  algorithm  creation  to  obtain  the  same  accuracy  for  the  minimal 
processing  time.  Allan  and  Barnes  in  [1]  showed  that  the  filtering  effect  strongly  depends  on  a  measuring 
time  interval  A  and  that  it  must  be  of  100...  1000  sec.  Later,  series  of  reports  have  been  devoted  to  the 
Kalman  filters  application,  especially  for  phase,  time,  and  frequency  real-time  estimates  with  prediction  [2- 
5].  The  papers  are  based  on  the  Kalman  fundamental  approach  [6]  for  the  discrete-time  optimal  linear 
signal  estimation  with  a  white  Gaussian  noise  and  developed  later  by  many  authors  in  [7-8],  for  instance. 

The  report  addresses  to  the  results  of  the  various  discrete-time  Kalman  algorithms  (one-,  two-,  and  three- 
state)  use  for  the  estimation  of  synchronization  errors  being  based  on  the  timing  signals  of  Motorola 
Oncore  UT+  Receiver.  Once  the  purpose  was  to  compare  the  various  Kalman  algorithms, then  the  results 
were  found  for  one  OCXO  unit  with  an  AT-cut  resonator.  While  doing  so,  we  dealt  with  many  measured 
and  estimated  functions  of  the  time  error,  freauency,  and  aging  related  to  the  OCXO  and  rubidium 
standard.  It  allowed  comparison  and  selection  of  appropriate  optimal  filter  structure  separately  for  each 
type  of  an  oscillator  and  for  measurement  ana  synchronization  tasks.  As  it  was  expected,  different 
algorithms  gave  different  filtering  errors.  Only  when  the  frequency  drift  showed  a  quasi-stationary  nature 
did  it  not  matter  what  type  of  Kalman  filter  was  used.  On  the  contrary,  in  the  non-stationary  case 
only  optimal  filtering  and  control  algorithms  provide  sufficient  accuracy.  Here  the  major  point  is  the 
timing  model  of  the  steered  source.  One  must  expect  that  a  rubidium  source  of  frequency  may  be 
disciplined  by  GPS  signals  with  relatively  small  error  in  comparison  with  the  crystal  one, because  of  good 
a  priori  prediction  of  the  phase  behavior  in  time.  We  study  major  of  these  possibilities  in  the  report. 
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DISCRETE-TIME  KALMAN  FILTER  MODEL 


The  oscillator  or  clock  discrete-time  equations  are  given  by  observation  and  state  vectors  [9] 

^v=HvXv+uv+n0v,  (1) 

=  Av-i^v-i +nxv’  @) 


where  v  =  0,  1,2,  ...  corresponds  to  discrete-time  tv  and  measuring  time  interval  A  =  tv  -  fv-i ,  ~  ^(N)  is 

m-dimensional  observation  vector  formed  by  the  reference  short-term  noisy  GPS  timing  signals  and  the 
oscillator,  Lv  =  \{tv)  is  ^-dimensional  oscillator  state  vector  (time,  phase,  frequency,  aging,  ...),  Hv  = 
H(/v)  is  mxn  dimensional  measurement  matrix,  uv  =  u(/v)  is  m-dimensional  vector  contains  the  control 
signals,  Av  =  A(rv)  is  nxn  dimensional  state  transition  matrix,  n0v  =  n0(/v)  and  nXv  =  nx(tv)  are  jointly 
independent  vector  white  noises  with  zero  expectations  and  covariance  matrixes  Vv  =  V(/v)  and  Tv  = 
'¥(tv)  are  of  mxm  and  nxn  dimensions  correspondingly 

Vv  =  £{n0Xv}>  (3) 

Vv  =  E{nXvni}.  (4) 


As  usual,  they  deal  with  single  observations  and  estimate  some  states  of  an  oscillator  or  clock.  It  means, 
by  expansion  that  if  m  <  n  then  one  may  use  the  following  algorithm  of  linear  Kalman  filtering  based  on 
(l)and(2) 


K  =  K-xK-i + Kv<£v  -  uv  -  hvav_,  v, )  ’ 

Rv  =(I-KVHV)RV, 

Rv  =  Av_,Rv_,Av-i  +  j 

Kv  =  RvHv(HvRvHv+Vv)-', 


(5) 

(6) 

(7) 

(8) 


where  X  =  X(tJ  is  a  vector  of  oscillator  state  estimates,  I  is  unit  matrix,  and  Rv  =  R(L)  is  the  error 

covariance  matrix.  Solution  (5)  is  justified  for  a  common  case  and  may  not  be  simplified,  as  a  rule. 
Nevertheless,  in  simple  particular  cases,  one-dimensional  algorithms  may  be  sufficient.  Setting  n  =  m  =  1 
in  (5) — (8)  we  come  to  the  following  general  form  for  the  estimate 


=  PvX,  +  *v(£v  -  «v  -  # vPv-l^v-1  )  > 

(9) 

i  i  ,  K 

(10) 

R, ,  Pv-A-.  +  Aw  A*’ 

k  =HV^~, 

(11) 

V  v  D 

•^Ov 

where  all  the  denotes  correspond  to  the  vector  case. 

ONE-DIMENSIONAL  ESTIMATES 


Based  on  (5) — (8)  we  consider  below  some  particular  algorithms  of  one-dimensional  Kalman  filters 
intended  for  frequency  estimation  only  in  crystal  and  rubidium  oscillators  if  GPS  timing  signals  are  used 
as  the  reference  ones.  All  results  given  are  obtained  with  the  Motorola  Oncore  UT+  Receiver. 

STATIONARY  KALMAN  FILTER 


Let  pv  =  p,  DXv  =  D-h,  Hv  -  H,uv-  u,  and  D0v  =  A>  be  the  constant  values  in  time  and  the  limit  Ru  =  fim 

v— ►  « 

be  known.  In  this  case,  the  following  equation  forms  the  Kalman  estimate 


(12) 
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which  is  the  stationary  and  non-optimal  one,  strictly  speaking.  Taking,  as  usual,  that  H  =  1 ,  (3  =  1 ,  and 
u  =  0  we  transfer  to  the  most  simple  form  of  the  stationary  Kalman  filter 

O3) 

where  R  is  defined  from  (10)  with  t  00  and  Rv  =  Rv_]  by  the  equation  R]  +  RvDXv  =  DXvD0v.  As 

follows  from  simplification,  such  filter  must  bring  small  error  in  the  steady  state  mode  only.  The  results 
obtained  by  the  filter  (13)  are  denoted  in  the  Figures  as  the  l"  estimate. 

Quasi-optimal  Kalman  filter 


Ifwetake  that  (3V  =  1  ,HV  =  1,  and  wv  =  0,then  we  come  to  the  following  solution 

K  —  ^v-i +  Kkv  —  ^v-i)  > 

1  _  1  1 

K  K-t +  £\jv 


(14) 

(15) 

(16) 


If  Z),v  =  Dx  and  D0v  =  D0  in  (14) — (16)  are  the  constant  magnitudes, then  the  quasi-optimal  Kalman  filter 
appears.  Here  we  may  take  the  following  assumptions.  Let  DXv  =  A.o  and  D0v  =  £>00  be  equal  to  the  values 
known  initially.  For  this  case,  we  obtain  the  function  marked  in  the  Figures  as  the  2nd  estimate.  We  also 
may  evaluate  the  variances  DXv  =  Dx  and  £>ov  =  Do,  based  on  the  regression  line  and  obtain  Kalman 
estimate  starting  from  the  certain  discrete  point.  In  Figure  1,  we  denoted  this  case  the  3rd  estimate. 


Optimal  Kalman  filter 


Optimal  estimates  are  obtained  by  (14) — (16)  if  all  the  necessary  parameters  of  the  observation  and 
frequency  behavior  are  known  for  the  arbitrary  time  point.  We  get  such  values  (strictly  speaking,  it  is  also 
quasi-optimal  ones)  based  on  the  regression  line 

y,  =  I  +— Ps,(',-0’  i  =  0’m  O7) 

where  ^  =  and  t  =  E{t,}  are  expectations  of  and  t:  respectively,  a]  =£{(£,,- ^)2}  and 

a]  =  E{{t ,  -/)2}  are  the  corresponding  variances,  p  =  is  the  correlation  coefficient. 

Procedure  (17)  takes  first  m  points  of  the  process  for  definition  of  the  filter  parameters  by  the  following 
way.  We  put  down  k(0)  =  kv_1  =  ym;then  we  find  out  d^v  =  -  yv  and  obtain  the  values  D0v  -  E{d£} , 

(18).  We  further  take  the  variance  of  the  regression  function  (17)  as  cr2  =  M{{yj  -y)2}  and  get 
DXv  =  +  o; ,  (19),  where  DXos  is  expected  variance  of  an  oscillator  frequency.  Let  us  note  that 

function  DXv  depends  on  the  regression  function  slope  that  results  in  the  filter  dynamic  properties  (inertia 
and  inaccuracy).  Finally,  we  take  the  initial  value  R^  {  =  Dyjn .  The  optimal  estimates  are  denoted  as  the  4,h 
ones  in  the  Figures. 

Experimental  Results 

Figures  1-3  show  the  observation  £v  formed  by  the  OCXO  frequency  averaged  over  the  time  interval  A 
=  100  sec  w'ith  the  GPS  timing  signals.  Here  we  also  show  the  measured  frequency  function  f,  obtained 
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with  use  of  the  rubidium  standard.  We  take  it  to  compare  with  the  received  Kalman  estimates  and  to 
evaluate  the  expected  filtering  errors  for  different  algorithms.  To  get  generalization  we  use  the  same 
observation  for  the  all  Kalman  filters,  noting  that  many  experimental  data  have  been  processed  in  reality. 
Both  l*'  and  2nd  estimates  are  based  on  the  same  constant  values  D0v  =  D0  =  5-10'19,  DXv  =  Df  =  10'22,  and 
start  zero  point  A,0  =  0.  In  opposite,  3rd  estimate  is  started  approximately  after  three  hours  of  regression 
processing  and  based  on  the  automatically  defined  values  X(0)s^v_,  =ym,  D0v  (18),  and  DXv  (19) 
obtained  by  the  linear  regression  function  (17).  Another  benefit  of  the  regression  function  (17)  usage 
jointly  with  the  one-dimensional  algorithm  (14) — (16)  is  the  possibility  of  optimal  Kalman  filter  creation. 
The  procedure  is  similar  to  the  above-considered  (17) — (19).  The  difference  is  following.  Here  we  use 
the  sliding  regression  estimates  of  Dqv,  DXv,  and  Lv-  based  on  which  estimates  A.v^  are  provided  through 

the  Kalman  algorithm  (14) — (16).  Figure  3  shows  the  optimal  function  (4*  estimate)  jointly  with  the 
curve  obtained  by  the  quasi-optimal  filter  (5*  estimate).  The  major  difference  between  those  two 
algorithms  is  following.  In  the  4*0356  we  use  the  variance  function  DXv  =  Dx  +  ZX  v ,  where  ZX  js  the  de¬ 
sired  filtering  effect,  and  Z\vis  the  additional  increment  depending  on  the  estimate  function  rate  so  that 
Au  =  (K  ~  A-iVA  ■  t^ie  case,  only  equality  DXv  =  Dx  is  substituted. 

TWO-DIMENSIONAL  ESTIMATES 


Let  us  address  to  the  two-dimensional  Kalman  filter  pursuing  the  aim  to  get  jointly  the  estimate  fv  of 


frequency  fv  and  estimate  &v  of  its  aging  rate  av,andto  obtain  the  prediction  fv.  Based  on(l) — (8)  we 
write  [9] 

4V  =  HAv  +nQv,  (18) 

*■»  =  Aw_iVi  +nxv,  (19) 


where 


Hv  =  [l  0], 


Av 


/v 

aj’ 

1  A' 
0  1 


(20) 

(21) 

(22) 


nQv  is  a  white  Gaussian  noise  with  discrete-time  covariance  V0v  =  D0v  =  £{«2V }  =  N0v/2A  and  nXv  is  a 


vector  white  Gaussian  noise  with  covariance  matrix  [8,  3]  qp  =  ^gv  /\ 

AV  2 


A2/3  All 
All  1 


,  where  Nnv/2  and 


Nav  /  2  are  two-side  continuous  spectral  densities  of  the  noises  n0v  and  «av .  Correspondent  prediction  is 


formed  through  estimates  /v  and  av  by  /v+1  =  /v+Aav.  Figure  3  shows  the  predicted  curve  (6'h 
estimate)  jointly  with  4,A  estimate,  for  comparison.  Here  we  see  that  the  error  in  transient  of  prediction  is 
more  appreciable  in  comparison  with  the  one-dimensional  optimal  case.  Nevertheless,  both  optimal 
estimates  obtain  rather  similar  errors  for  the  stationary  observation  data  through  20  hours.  So,  in 
principle,  direct  frequency  estimate  is  sufficient  for  prediction;  that  is  extremely  important  for  the 
synchronization  needs  when  the  GPS  timing  signal  is  not  available. 


THREE-DIMENSIONAL  ESTIMATES 


Here  we  address  to  the  most  accurate  Kalman  algorithm  based  on  the  oscillator  timing  model  [10]. 
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Mathematical  model  of  the  oscillator  timing  signal 


In  the  ideal  case  the  total  phase  <X>id(7)  of  an  ideal  timing  signal  is  presented  by 

a>ld(/)=27rvnomr  (23) 

where:  vnom  is  called  nominal  frequency.  Correspondingly,  total  instantaneous  phase  model  of  actual 
timing  signals  <&(t)  is  modeled  as: 

Q(?)  =  <D0  +  2nvnom  (l  +  y0)t  +  %DvnoJ2  +  cp(t)  (24) 

where:  d>0  is  the  initial  phase  offset,  y0  is  the  fractional  frequency  offset  from  the  nominal  value  vnom 
(mainly  due  to  finite  frequency  settability  of  the  clock);  D  is  the  linear  fractional  frequency  drift  rate 
(basically  representing  oscillator  aging  effects);  cp(r)  is  the  random  phase  deviation  component.  Based  on 
the  definition  of  time  error  and  the  above  model  (24),  the  following  model  for  time  error  x(t)  results: 


D-D. 


ref  ?2  + 


2nv„„m 


(25) 


Assuming  that  for  the  measurement  of  x(t)  the  independent  clock  configuration  applies  and  that  the 
reference  clock  is  properly  chosen  (i.e.  all  its  degradation  sources  yo.ref ,  Aef  >  and  cpref(0  are  negligible  as 
compared  to  those  of  the  clock  under  test),  the  x(t)  model  reduces  to: 

< v(t ) 


,  ,  D  2 

-  xo+yo{  +  n  1  +  o 

2  2;rv„, 


(26) 


When  the  synchronized  clock  configuration  applies  and  all  slave  clocks  involved  in  the  distribution  of 
timing  (including  the  clock  under  test)  are  operating  in  locked  mode,  y0>ref  =  yo  and  Aef  =  D  can  be 
assumed;  the  x(t)  model  then  reduces  to: 


x{t)  =  x0  + 


<PU)  - 


2  7tvn 


(27) 


Three  dimensional  Kalman  filter 


Based  on  the  oscillator  timing  model  (26)  and  the  corresponding  Kalman  filter,  one  may  expect  the  best 
filtering  effect  if  the  model  and  the  filter  are  matched.  Let  us  consider  this  possibility  in  detail  using  the 
model  (26)  and  the  algorithm  (5) — (8)  to  create  the  corresponding  three-dimensional  Kalman  filter. 
Decomposing  the  time  error  function  xv  into  the  series,  we  write 

A.  =  *v-i  +  3\.-ia  +  0.5av_,A2  +  ,  (28) 

where  yv  =  yv_,  +  av_,A  +  is  a  frequency,  and  av  =  av_,  +  nav  is  a  linear  component  of  frequency  aging. 
It  allows  us  to  write  the  equation  of  the  filter  state  in  a  form  of  (2)  where 


xv 

’l 

A 

A2  /2~ 

y. 

’  4  = 

0 

1 

A 

’  nxv  = 

_av_ 

0 

0 

1 

_ 

_  O.V  __ 

(29) 


One-dimensional  observation  equation  (1)  is  defined  as 

!U  =  fev]>^  =  [l  0  0],  n0v  =[n0v]  •  (30) 

Correspondingly,  the  Kalman  filter  equations  are  defined  by  (5) — (8),  where  for  independent  and 
incorrelated  noises  we  get 

K  =  E{noy0v}  =  D0v  =  SQv  i  =  ^  =  cr2ov ,  (31) 

A  2A 
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xllv 

vp 

x12v 

x13v 

=  £KXv}  = 

1  2iv 

vp 

1  22  v 

^23v 

? 

(32) 

L  X3lv 

vp 

x32v 

^33vJ 

rAe 

/3  +  SavA5 

/20;  ¥21, 

=  VT/ 
T12v 

< 

n 

1 12  + 

OO 

5 

co5 

T22v  =  ^A  +  5a»A3/  3; 

TV  : 

J  JV 

=  s*. A ;  %u  =  ^  = 

SavA3/6; 

SavA2/2; 

V2,  = 

Ny./  2  and  S„  =  NaJ  2 

are  correspondent  two-side  spectral  densities  of  the  continuous  noises  depending  on  the  time  interval  A  in 
a  general  case.  If  to  account  an  aging  noise,  only  then  matrix  (32)  yields 


N 

¥  =f-s^A 
v  2 


‘A4/ 20 

A3/8 

A2/ 6 

=  ^A 

'0 

0 

o' 

A3/ 8 

A2/ 3 

A/2 

0 

0 

0 

2 

A2/ 6 

A/2 

1 

0 

0 

1 

(33) 


Corresponding  prediction  equations  are  obtained  by  the  above-considered  algorithm  being  based  on  the 
estimates  of  xv ,  yv  and  av  as  follows 

=  xv  +  yv  A  +  0.5avA2 , 

L,=jv  +  «v  A,  (34) 

av+1  =  av . 

The  equations  (34),  at  bottom,  are  the  working  ones  for  the  holdover  operation  mode  of  the  slaved 
oscillator  when  GPS  signals  are  not  available  [11], 


Experimental  results 


The  results  are  given  below  for  the  Kalman  algorithm  (28) — (34)  use  for  the  synchronization  error 
estimation.  We  consider  the  OCXO  phase  (time  error)  measured  for  A  =  1 00  sec  based  on  the  GPS  timing 
signals  of  the  Motorola  ONCORE  UT+  Receiver.  Here  digital  sequence  xv  (the  oscillator  phase  reduced 
to  the  time  interval  error)  is  the  initial  time  error  data  that  is  processed  by  three-dimensional  Kalman  filter 
(28) — (34).  As  a  result,  we  get  the  estimates  xv,  yv,  and  av.  Frequency  estimate  yv  is  compared  with  its 

accurate  magnitude  obtained  by  direct  measurement  with  the  use  of  the  reference  time  signal  of  lsec  of 
the  rubidium  standard  in  the  same  time  scale.  Figure  5  shows  the  initial  observation  (time  error)  xv  jointly 

with  its  Kalman  estimate  xv  and  expected  error  calculated  as  the  difference  s„  =  xv  -  xv .  Flowing  from 
the  physical  sense,  spectral  density  of  the  fluctuation  of  the  aging  speed  within  the  noise  matrix  (33)  is 
taken  as  Sa0  =0.5Aa0A  =  4-10'‘2Hz2  for  the  process.  With  this,  the  expected  root-mean-squire  time 
deviation  (TDEV)  and  maximal  time  interval  error  (MTIE)  are  estimated  as  cE  =  282.7  ns  and 
Axmax  =  39680  ns.  Thus,  in  principle  we  get  the  possibility  to  reduce  the  time  error  by  digital  PLL  to 
Ax  /a  =140  times.  Figure  6  exhibits  the  measured  function  of  the  relative  frequency  behavior  of  the 
same  OCXO  jointly  with  the  Kalman  estimate.  Ifwetake  the  measured  data  as  the  reference  data,  then  the 
difference  function  may  be  taken  as  the  error  function,  where  TDEV  equals  oEr  =  146-1  O'17  for  the 

considered  case.  Correspondingly,  TDEVs  of  the  measured  and  estimated  curves  are  a  =  325  •  10"12  and 
av  =  286  •  1 0~12 .  Hence,  the  expected  frequency  control  effect  yields  a,  =  a  /a  =  1 .96  and 
avm  =  /oty  =  2.4  with  respect  to  the  measured  and  estimated  data  correspondingly.  Finally,  Figure  7 
shows  estimates  obtained  for  the  frequency  aging.  Mean  value  of  av  equals  av  =  £{a(}  =  -2.7  •  1CT15  /s 
upon  the  observation  interval.  TDEV  for  av  is  equal  to  oa  =  3.4  •  10~I4/c,  and  the  ratio  of  the  two 
statistical  estimates  is  tfa  /av  =  12.4  that  is  due  to  the  extremely  small  SNR«1. 
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CONCLUSION 


We  have  considered  several  Kalman  algorithms  intended  for  the  time  error,  frequency,  and  aging 
estimation  of  the  slaved  oscillator  being  based  on  the  Motorola  GPS  Oncore  UT+  Timing  Receiver. 
General  findings  of  the  studies  are  following.  The  best  results  of  the  parameter  estimation  are  obtained  if 
the  timing  model  of  the  oscillator  is  known  and  matched  with  the  filter  dimension.  Thus,  selection  of  the 
proper  Kalman  algorithm  has  to  be  done  at  the  early  stage  of  the  filtering  and  control  to  get  the  best 
effect.  In  this  regard,  three-dimensional  algorithm  is  the  most  efficient  for  the  rubidium  oscillator  and 
brings  less  good  results  for  the  crystal  one,  as  the  function  of  the  crystal  oscillator  frequency  behavior  in 
time  has  the  complex  nature  and  is  less  predictable. 

It  follows  that  the  stationary  one-dimensional  filter  obtains  accurate  estimate  with  rather  great  noise.  The 
quasi-optimal  one  with  unknown  initial  statistics  shows  big  error  at  the  beginning  stage  of  the  process. 
Approximately  the  same  results  give  the  optimal  one-  and  two-dimensional  filters.  Therefore,  one¬ 
dimensional  filters  seem  to  be  preferable  for  the.  fast  and  relatively  accurate  estimate  of  the  only 
parameter  (frequency,  for  instance).  Two-dimensional  ones  allow  prediction  of  the  freauency  offset  and, 
therefore, are  more  preferable  for  synchronization  tasks.  Three-dimensional  Kalman  filters  exhibit  small 
inertia  (negligible  dynamic  error)  and  hieh  accuracy  jointly  with  prediction  of  the  time  error  for  the 
holdover  operation  mode.  This  filter  is  well  matched  with  a  time  model  of  the  rubidium  standard,  and 
increase  of  the  filter  dimension  is  required  in  the  case  of  a  crystal  oscillator,  as  usual.  The  filter 
adaptation  gives  good  results  for  one-  and  two-dimensional  algorithms.  Contrarily,  high  .adaptation 
efficiency  for  three-dimensional  filters  take  place  in  the  case  only  when  the  time  error  function  exhibits  a 
step  character.  Let  us  mark  that  such  a  case  is  not  the  usual  one  for  the  rubidium  oscillators. 
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Figure  1.  The  Kalman  estimates  of  the  OCXO  frequency,  GPS-based  observation,  and  measured  data. 


Figure  2.  Magnified  scales  correspond  to  estimates  obtained  by  the  Figure  1 . 
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Figure  3.  Magnified  scales  correspond  to  4th  and  5lh  estimates  and  measured  data. 
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Figure  4.  Magnified  scales  correspond  to  4th  one-state  estimate  and  6“  two-state  estimate  (prediction). 
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Figure  5.  GPS-based  three-state  Kalman  estimates  of  time  error  xev  =  xv  (a),  frequency  offset  yv  =  y 
and  measured  data  (b),  and  frequency  aging  rate  av  =  ctv  (c)_ 


